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Abstract 

We present an implementaion of interface between the full-potential linearized augmented plane wave package 
Wien2k and the wannier90 code for the construction of maximally localized Wannier functions. The FORTRAN 
code and a documentation is made available and results are discussed for SrVC>3, S^IrCU (including spin-orbit 
coupling), LaFeAsO, and FeSb 2 . 
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1. Introduction 

In recent years, we have seen a revival of Wannier's idea [1] to express the single electron excitations 
in a periodic potential in terms of "orthogonal 'atomic' wave functions" , later coined Wannier functions. 
While the calculation of bandstructures as, e.g., within density functional theory (DFT) [2] is most efficiently 
performed in the basis of extended wave functions indexed by reciprocal-space vectors k, the complementary 
Wannier picture is often useful - if not essential. The latter does not only provide insight into the nature 
of chemical bonding, but the localized Wannier orbitals are more suitable to describe physical phenomena 
where the local aspect is important. These range from the dielectric polarization [3], e.g., in the context 
of ferroelectrics [4,5], and orbital polarization [6] to molecular dynamics [7] and transport through nano 
structures [8,9]. 

Also to describe electronic correlations beyond the local density approximation (LDA) to DFT, Wannier- 
like orbitals are very helpful. Here they allow us to limit the interactions to a computationally manageable 
subspace, such as the local interactions between electrons on the same transition-metal atom in the d- 
orbitals. The on-site term represents not only the largest constribution to the interaction, but typically gives 
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raise to the dominant electronic correlations. This is because inter-site interaction terms yield mainly the 
Hartree contribution if the number of neighbors is large [10], and this is already included in LDA. The 
concept of Wannier functions has hence found its use also in materials specific many-body calculation with 
the LDA+dynamical mean field theory (DMFT) approach [11]. The first LDA+DMFT calculations made 
use of downfoldcd linearized muffin tin orbitals (LMTO) [12] and Nth order muffin tin orbitals (NMTO) 
[13] which exploit the local character of the basis functions. For other bandstructure methods with bases of 
extended orbitals a transformation to Wannier functions is needed. 

While the conversion from Bloch to Wannier functions is essentially a Fourier transformation, there 
remains an arbitrariness in choice of the phases of the Bloch functions; i.e., this conversion is not unique. 
This arbitrariness can be used to construct localized orbitals. Marzari and Vandcrbilt [14,15] developed a 
procedure to fix the phases (and further k-dcpendcnt unitary transformations) so that the resulting Wannier 
orbitals minimal spread around their centers. Together with coworkers, they developed a corresponding 
program package wannier90 [16]. 

Since then the maximally localized Wannier functions (MLWF) have been implemented for many DFT/LDA 
program packages ranging from plane waves to LMTO [17,18] and linearized augmented plane waves (LAPW) 
in the FLAPW [19] and FLEUR code [20]. Alternatives to MLWF have been proposed as well, including the 
recent Wannier function projection [21,22] where the overlap of the Bloch functions with the d- or /-orbitals 
within the atomic or muffin tin spheres defines the Wannier functions. This approach has been implemented 
for LMTO [22] and by Aichhorn et al. [23] for the FLAPW of the Wien2K code [24]. While this approach 
is somewhat easier to implement, the final Wannier functions are no longer uniquely defined and may differ 
considerably. With a clever choice of the functions to project onto, results very similar to MLWF can be 
obtained. 

In this paper, we present an interface [25] between Wien2K [24] implementation of the FLAPW method 
and the wannier90 code [16]. This iterface generates the wannier90 input files (.eig, .mmn, .ami) from the 
Wien2K electronic structure. Futhermore, we provide a code for the generation of MLWFs in direct space, 
primarily useful for visualization. The interface consists of FORTRAN90 programs (tested for gfortran 
and ifort; hitherto not parallelized) with additional shell and python scripts to simplify the work-flow. A 
manual is available at [25]; the work-flow will be briefly described in Section 2.4 The runtime of the main 
program (w2w) was 4.6 s for a simple test system such as SrV03 with 64 k-points, running under Linux 
with a 64-bit Intel Xcon 3 GHz 4GB RAM. For visualization, an output to XCrysSDen [26] is also provided. 

The following sections are organized as following: In Section 2, we discuss the theoretical background and 
equations, in particular the projection onto maximally localized Wannier function (MLWF) in Section 2.1 
and the LAPW of Wien2K in Section 2.2. Some practical considerations are presented in Section 2.3. The 
workflow of our Wien2Wannier interface is outlined in Section 2.4. In Section 3, we present exemplarily some 
applications of the code. We start with the simple bandstructure of SrVOa in Section 3.1 to demonstrate how 
the choice of the energy window affects the shape of WFs and the corresponding tight-binding Hamiltonian. 
Sr 2 Ir04 with an idealized structure illustrates in Section 3.2 the construction of WFs in a system with strong 
spin-orbit coupling leading to orbitals with complex phases. In Section 3.3, LaFeAsO is used to show the 
utility of WFs for backfolding of the electronic bandstructures, in order to obtained the simplest possible 
k-space representation of the electronic structures. Finally, in Section 3.4, we show the application of the 
code to a material with low site-symmetry by the example of FeSb 2 . 

2. Theory 

2.1. Wannier functions 

Following Refs. [14,15] we define the Wannier orbitals u> n R,(r) in the i?'th unit cell through a unitary 
transformation of the Bloch eigenstates tp n k( r ) of the periodic Hamiltonian: 




(1) 
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In the following the overall phases of the Bloch functions are assumed to be arbitrary, except of the periodic 
gauge constraint VVnk+G = tp m k for G being a reciprocal lattice vector. According to the Bloch theorem the 
function ^nk factorizes into a product of the plane wave exp(zk • r) and a periodic function u^: 

V-nk = e lk ' r u nk . (2) 

We are looking for that £/( k ) which yields MLWFs, i.e., minimizes the spread (variance) 

^2(w mR (r)\r 2 \w mR (r)) - (w m R.(r)|r|w mR (r)) 2 . (3) 

m 

As the Wannier orbitals are expressed in terms of the Bloch wave functions in Eq. (1), we need to calculate 
the following overlap integrals on a uniform k-mcsh in the Brillouin zone (BZ): 

= (u mk \u nk+b ) (4) 

= (^ m k|e^ b ' r |^„k + b> (5) 

For further details on how to determine the MLWF from these overlap integrals, we refer the reader to 
[14,15]. 

Moreover, for the disentanglement procedure described in [15] and a good starting guess, also overlap 
integrals with a set of hxed trial orbitals g n 

A%1 = (Vmkk) (6) 

have to be calculated. The evaluation of and Ami in the LAPW basis is the subject of the rest of 

this section. 

2.2. Linearized augmented plane-waves 

In the LAPW method, the direct space is divided into two regions: the interstitial space I and the non- 
overlapping muffin-tin (MT) spheres around the nuclei at R^. The Bloch functions are expanded into 
plane waves in I and into partial atomic waves in the MT spheres 

Ij^G)^ pel 

K(*)=< L, G a _ B . a . „. „. B . x (7) 
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where Yg m are the spherical harmonics and u±, Ui, and 112 are numerically determined functions of the radial 
variable rp — |r — R^|, which are defined separately for each MT sphere (3, orbital quantum number I 
and, in spin-polarized calculations, spin projection a. For simplicity wc do not show the dependence of the 
expansion coefficients C£, a\^ im , b^ tm , and c\^ lm on the band index. Both LAPWs and augmented plane 
waves plus local orbitals (APW+lo) basis sets can be used for our interface. While further details can be 
found in Ref. [27] we point out two features of the expansion: (i) it is not possible to obtain the periodic 
part it n k inside the spheres in a simple straight-forward way from the above expansion of the Bloch function 
and therefore we used Eq. (5) in actual calculations; (ii) the Bloch function is uniquely defined by C£(G) 
since the other expansion coefficients are in fact functions of C£(G), defined by the continuity conditions at 
the surfaces of the MT spheres. 

The calculation of Mmn°^ naturally splits into the interstitial part and the MT part. To evaluate (5) inside 
the MT sphere /?, the plane-wave exp(— ih ■ r) is expanded into products of spherical harmonics and Bessel 
functions using the well known formula 

e- ib - r = 4^ e - jb - R ^ i ^(6r / ,)^f, m (-b)y £m (f ). (8) 

l m 

Substituting Eqs. (7) and (8) into Eq. (5), one obtains an expression for the contribution of within 
the MT spheres in terms of products of radial integrals J^ 13 dr r 2 Up Cr ei (r)j£ 2 (br)u^' 7 i3 (r); angular integrals 
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((Yt 1 miYe 2 m2Ye 3m3 )) , known as Gaunt numbers; and an additional weight factor given by the corresponding 
products of the expansion coefficients. 

The interstitial contribution Mmn^ on the other hand reads 

^(GO^VbCG)^^- ')"-),, (9) 

GG' 

where the ()/ integration runs over the unit cell with the MT spheres removed, which gives 
{e HO-o').r h = v6g g , Wfj J- e -HO- G 'yn p sin(x) - , cos(x) 

p 4 x (10) 

where x = pp\G - G'|, Vp = -7rp|, 

and pp is the radius of the MT sphere at Hp. 

The periodic gauge has to be enforced when the integration in Eq. (1) wraps around the Brillouin zone 
boundary, namely a reciprocal vector G must be chosen such that k + b G lies within the Brillouin zone. 

The construction of Ami, is rather simple if we allow trial orbitals to be non-zero only inside the MT 
sphere, i.e. of the form 

9n =YY J 9i m {n)u^ l {rp)Y lm {rp), (11) 

p Cm 

where the list of non-zero coefficients g^ m (n) is provided as input. Using (7) we obtain 

= E E (<a>?>?a + 6& m <«?>?A + <^<</i</>r) (12) 

p Im 

where () r denotes the radial integral with r 2 dr. 



2.3. Practical considerations 



In the following we briefly describe the steps taken in computing WFs. We start by generating a uniform 
k-mesh in the Brillouin zone (running KGEN with only the identity matrix in the list of symmetry opera- 
tions) followed by generating the eigenstates using LAPW1. In case of calculations with spin-orbit coupling 
LAPWSO should be executed. Spin-orbit calculation must be run as formally spin-polarized even if the 
polarization is zero. Before running the wien2wannier interface code the initialization run of wannicr90 must 
be executed to generate the case.nnkp file. The inputs to wien2wannier specifies the bands in the initial 
Hilbert space, the number of the bands in the target Hilbert space and the expansions (11) of the target 
orbitals. The wien2wannier generates the necessary files case.eig, case.mmn and case.amn, which serve as 
an input for the MLWF construction by wannier90. In case of spin-orbit coupled calculations, a separate 
run of wicn2wannier is executed for each of the two components of the spinor wave function. The corre- 
sponding elements of the resulting case.mmn(amn) files are added up to form the total overlap matrices 
between the spinor functions. The postprocessing of the wannier90 results is, to a large extent, independent 
of the bandstructure code. However, for plotting the WFs in direct space information on underlaying basis 
functions is necessary. We have modified the LAPW7 code of the Wien2K package to generate direct space 
WF maps, which uses the transformation between the initial Bloch states and resulting WFs extracted from 
the wannier90 calculation. 



2.4. Computational details 

The calculations reported in the following Section were performed with Wicn2k [24] code employing the 
LDA exchange-correlation functional for SrVC>3, Sr 2 Ir04, LaFeAsO and the generalized gradient correction 
(GGA) for FeSb 2 ; in all cases the APW+lo basis set [28] was employed. To construct the Wannier orbitals 
we followed the sequence: (i) self-consistent bandstructure calculation using the irreducible part of the BZ; 
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(ii) generation of a uniform k-mesh throughout the entire BZ (wien2k) and of the corresponding k, k + b-list 
(wannier90); (iii) calculation of the Bloch eigenstates and eigenvalues on the new k-mesh (wien2k); (iv) 
evaluation of the M^l b ' and Ami elements (wien2wannier); (v) generation of MLWFs (wannier90); (vi) 
post-processing of MLWFs to generate tight-binding bandstructures, hopping integrals, plots of Wannier 
orbitals (wannier90+wien2k). In case of calculations with spin-orbit coupling the Bloch functions have two 
non-zero components, indexed with the spin quantum number, with a definite mutual phase. Once the two 
components are generated, the contributions to Mm^ and Ami of each component arc computed separately 
along the lines of the previous section and added up at the end. To plot the Wannier orbitals we extract 
the Uml matrices of Eq. (1) from the wannier90 processing and combine it with the wien2k utility for the 
generation of the Bloch state on direct space grids. 



3. Results 



3.1. SrV0 3 



SrV03 is a thoroughly studied material with a rather simple band structure (see Fig. 1) consisting of 
isolated groups of bands derived from O-p, V-d-t 2g and V-d-e g orbitals. Therefore it is very well suited as a 
testing case. We focus on the V-d-t 2g states, which (being partially filled) are of most interest. We compare 
two different settings, i.e., consisting of two choices of the Hilbert space to be represented by the Wannier 
functions. Namely, the space spanned by (i) t 2g bands only and (ii) all the V-d and O-p bands. In both 
cases, the t 2g Wannier functions are orthogonal to the O-p Bloch or Wannier states. In the former case (i) 
also the hopping integral between the t 2g and O states is still zero at the price of more extended Wannier 
orbitals with a substantial weight at the O sites. Naturally, the larger energy window (ii) clearly allows the 
construction of more localized WFs as demonstrated in Fig. 2. The Wannier orbitals (i) have visible density 
on the neighboring O sites, reflecting the mixed character of the t 2g bands. In contrast, the orbitals (ii) 
constructed from the larger energy window do not have any appreciable (visible) density at O sites (on the 
scale of the figure) since the latter would be now assigned to the O-p orbitals, which are explicitly presented 
by O-centered Wannier orbitals. 

The spatial extent of the Wannier orbitals is also reflected in the hopping integrals, which were calculated 
by the wannier90 code. With the small energy window (i) the t 2g bands are well described when at least nn- 
and nnn-hoppings is considered (see Fig. 1), which in the {xy, yz, zx} basis read 

/-268 00\ / 7 10 ^ 



tioo[meV] = 









-30 



tioi[meV] 



Hi 7 ii 



(13) 
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The remaining directions follow from symmetry considerations. The longest nnn t 2g -t 2g hopping corresponds 
to a length of of 5.4 A. Using the more localized orbitals (ii) we can achieve similar accuracy (see Figure 1) 
by considering only V-V nn-hopping and V-0 nnn-hopping, which translate into a direct spatial cut-off of 
only 4.3 A. The obvious price to be paid are larger matrices (14 x 14 for (ii) vs 3 x 3 in case (i)). We point 
out when no spatial cut-offs are introduced, both choices (i) and (ii) represent the t 2g bands to the same 
arbitrarily high accuracy, determined by the size of the k-mesh used to construct the MLWFs. 



3.2. Sr 2 Ir0 4 



Sr 2 Ir04 has been recently the subject of intense investigations due to the close connection between the 
spin-orbit coupling and its Mott-insulating ground state. We use it as a example of a material where spin- 
orbit coupling substantially modifies the band structure and leads to Wannier orbitals, in which both spin 
projections are mixed. For sake of simplicity we have performed the calculations using an idealized double- 
perovskite structure, while the real material is characterized by a tilting of the IrC>6 octahedra. In the 
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Fig. 1. (color online) LDA band structure of SrV03 (circles) with dominant contributions marked with color: O-p (blue), 
V-d-t2 g (red), and V-d-e 3 (black). Tight-binding bands obtained from MLWFs for (i) [left panel] only the t2 g band (ii) [right 
panel] all V-d and O-p bands, with spatial cut-offs as described in the text (the Fermi level is at 8.16 eV). 




Fig. 2. (color online) LDA band structure of Sr2ir04 (back circles) together with the MLWF tight-binding fit to the J = 1/2 
bands for different direct space cut-offs; the Fermi level is at 8.2 eV. (left) The xy Wannier orbital plotted as an isosurface of 
the charge density |«)(r)| 2 and colored by the sign of ui(r)2. The left panel corresponds to the large energy window (ii) and the 
right panel to the small energy window (i) for the same isovalue. (right) 
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Fig. 3. (color online) LDA band structure of Sr2lr04 (back circles) together with the MLWF tight-binding fit to the J = 1/2 
bands for different direct space cut-offs (the Fermi level is at 8.2 eV). 

following, we discuss the LDA band structure of Sr2lr04, which represents a crude approximation to cor- 
relation effects due to local interactions. Within this picture the electronic structure (see Fig. 4) can be 
understood by considering the crystal-field splitting, spin-orbit coupling and the inter-site hopping via the 
manifold of Ir-e? bands. The crystal-field splitting, being the largest of the three, opens a gap between the 
tig and e g bands, rendering the latter empty, while the formed accommodate one hole per Ir atom. The ti g 
orbitals may be labeled with a pseudo-spin 7 = 1. The spin-orbit coupling further splits the ti g manifold 
into a quadruplet and doublet with pseudo-spin J = 3/2 and J = 1/2, respectively. Since the spin-orbit 
splitting is rather large the inter-site hopping leads only to moderate mixing of the states with different J. 
Therefore we may expect the isolated band doublet at the top of the t 2g manifold to be predominantly of 
,/ = 1/2 character. We hence construct the MLWFs for these two bands. 

For testing purposes, we have prepared Bloch states on k-grids of sizes 3x3x3,5x5x5, and 7x7x7. 
As trial functions we have used the J e ff = 1/2 orbitals, which adopt the following spinor form with respect 
to the local coordinate axes pointing towards the O atoms: 





10+) ~ > l*-> ~ > ( 14 ) 



where we omit the unimportant normalization factor. The results for the three different k-grids were almost 
identical with minor deviations only for the smallest grid. Due to a clever choice of trial functions, the original 
spread of the Wannier orbitals as defined in wannier90 actually changed by less than 0.5% during the MLWF 
optimization. This is consistent with the fact that the resulting MLWFs are essentially J e ft — 1/2 functions 
as shown in Fig. 4. In particular, the Wannier orbital, to a very good approximation, consists of a real xy 
orbital in one spin channel and a complex (x ± iy)z orbital in the other (±) spin channel. Note that the 
relative phase of the two components is not arbitrary and the corresponding charge density (sum of the spin 
components) has (approximate) cubic symmetry as expected for a J e ff = 1/2 orbital. 

In Fig. 4 we show convergence of the tight-binding bands to the original band structure with increasing 
maximum hopping distance. The band dispersion is governed by the nearest- and next-nearest-neighbor 
hopping and the bands are essentially converged when five coordination spheres are considered, which 
includes also the out-of-plane hoppings giving rise to the z-axis dispersion. The hopping amplitudes are 
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Fig. 4. (color online) Ir J = 1/2 Wannier orbital w+,[o,o,0] ( r ) visualized as a \w\ 2 isosurface. The almost real 4-spin component 
is blue, the "("-spin component is colored with cosine of its phase (red=real positive, green=real negative, yellow=imaginary), 
see color legend bar. 
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Fig. 5. (color online) The tight-binding hopping integrals of the one-band model as a function of inter-site distance for different 
BZ grid sizes (the lines are guides for eyes). 
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Fig. 6. (color online) LDA band structure of LaFeAsO (red cross symbols) compared with tight-binding bands obtained from 
MLWFs (green solid lines) for the dpp model (left) and the d model (right); the Fermi level is set to zero. 

summarized in Fig. 5. 
3.3. LaFeAsO 

In 2008, Kamihara et al. discovered superconductivity in F-doped LaFeAsO [29], immediately followed 
by the same finding in related materials. This opened a new field of research inviting many theoretical 
approaches to be applied, which in turn called for simplified models of the electronic structure capturing 
the essential chemistry of these systems. These new superconductors commonly have two-dimensional iron 
pnictide or iron chalcogenite layers, and the low-energy electronic structure around the Fermi level (Ep) 
consists of heavily entangled bands dominated by the Fe 3d states (see Fig. 6). While the multi-band character 
renders tight-binding fit impractical, the exact transformation using Wannier functions provides a controlled 
way of constructing tight-binding models representing the bandstructure to arbitrary accuracy, see Fig. 6 and 
[30-37]. 

Two types of tight-binding models have been considered for iron pnictides: the d model on the basis of 
Fe 3c?-like orbitals, and dp or dpp model using the Fe 3<i, pnictogen/chalcogen p, and oxygen p orbitals. In 
both cases the corresponding bands are reproduced exactly, as shown in Fig. 6, but the behavior of the two 
models differs if electron-electron interactions are explicitly included. As in the case of SrV03, the shape 
and spread of the Fe 3d orbitals depends on whether the d or dp model is used, see Fig. 7 and Table 1. 

One of the advantages of a tight-binding model, which represents the band dispersion of ab initio cal- 
culation, is that we can sometimes "unfold" the Brillouin zone. The unit cell of LaFeAsO contains two Fe 
and As atoms, where one of the As atoms sits above and the other below the Fe plane. Therefore, there are 
ten Fe 3d bands. It turns out that the translation symmetry of the tight-binding model is higher and the 
bandstructure can be unfolded to a five band scheme, corresponding to an effective unit cell with one Fe site. 
The corresponding bandstructure is shown in Fig. 8. The missing translation symmetry between the two Fe 
sites of the crystallographic unit cell means that the WFs on the two sites are not connected by a simple 
translation, but by a more general symmetry operation. Obviously, not all operators have the property of 
the tight-binding Hamiltonian that this does not matter; and therefore the five band model cannot replace 
the ten band one in general. 
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Table 1 

Spread of maximally localized Wannier functions for d model and dpp model (in units of A). 




Fig. 7. (color online) Isosurface contours of the x 2 -y 2 MLWO for the dpp model (left) and the d-only model (right). The 
amplitudes of the contour surface are +0.07 A~ 3 / 2 (red) and -0.07 A~ 3 / 2 (blue). 

In table 2, we list the hopping parameters of the d model. They are very similar to previous studies such 
as Ref. [30], Ref. [35], of Ref. [37]. It should be noted that since Fe 3d bands hybridize with La 5d, there can 
be slight difference, depending on whether we use the so-called "inner" window or not, see [37]. 
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Table 2 

Hopping integrals t(Ax, Ay; fi, v) in units of mcV. [Ax, Ay] denotes the in-plain hopping vector being different for different 
columns and (/t, v) the orbitals being different for the rows. The last three columns a y , I, and cr^ denote the symmetry 
transformations necessary to calculate t(Ax, — Ay; fi, u), t(— Ax, — Ay; fi, v), and t(Ay, Ax; /i, u), respectively. Here '±' means 
that the corresponding hopping is equal to ±t(Ax, Ay; fi, i/) in the same row, and '±(//,i/)' states that the hopping equals 
±t(Ax, Ay; ix 1 , u') in another (fi' , i/) row. Note that there is another symmetry relation t(Ax, Ay; fi, v) = t(—Ax, —Ay; v, p.). 
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Fig. 8. Band dispersion of the dpp (left) and d model (right) for LaFeAsO with the first Brillouin zone being extended to the 
second Brillouin zone, so that the unit cell contains only one Fe atom. Both Y and T' points correspond to T in the original 
band structure in Fig. 6. 

3.4. FeSb 2 

FeSb2 is a small band gap semiconductor which is intensively studied [38] for its unusually good ther- 
moelectric properties. In order to describe the Coulomb interaction within the Fe-3d shell and the strong 
hybridization between the Fe-3d and the Sb-4p electrons on an equal footing, a Wannier projection onto a 22 
band pd subspace is required. FeSb2 has an orthorhombic structure belonging to the space group 58(Pnnm) , 
in which the Fe sites are located at the centers of largely distorted Sb octahedra. This distortion results in 
a low point group symmetry at the Fe site. 

The LDA bandstructure, shown (grey solid lines in Fig. 10), is rather complicated, with about 10 bands 
in the energy interval of ±2eV around the Fermi level. The Wannier functions were constructed from 41 
Bloch states on a 10x10x10 k— grid, by projection onto the subspace of 22 MLWF in the energy region 
[-6eV,+6eV], using the disentanglement procedure of wannier90. 

We performed Wannier projections with different frozen energy windows and found that already the 
window [-6eV,3.9eV] allows for a good description of all the 22 bands in the whole energy interval of ±6eV 
around the Fermi level. The corresponding MLWFs are well localized with the spread ranging from 0.8A 2 to 
1.2A 2 for Fe-3d orbitals and from 3.6A 2 to 5. OA 2 for the Sb-4p orbitals. The largest spread of ~ I.2A 2 was 
obtained for the two Fe-e g -like Wannier functions which point more towards the direction of the ligands (see, 
e.g., first panel in the first row of Fig. 9), while the other three (Fe-^g-hke) Wannier orbitals are slightly 
more localized with a spread of 0.8A 2 -j- l.OA 2 (second and third panel of Fig. 9). The pictures of the three 
Sb-Ap orbitals are reported with the same scale in the second row of Fig. 9. The contributions of the different 
MLWFs to the bandstructure is shown as a so-called fat-band plot in Fig. 10, revealing a strongly mixed 
d-p character of all bands. 

By truncating the Hamiltonian in the MLWF basis set in direct space, we find that retaining hoppings 
up to 7.3A provides an accurate description of the LDA bands. A closer analysis shows that the strongest 
hopping processes (1.2 -f- 1.8 cV) are in the p-p sector of the Hamiltonian. Relatively large hoppings are 
found also in the p-d sector (up to 0.6 -j- 0.9 eV) while the d — d hopping amplitudes reach the maximum 
value of about 300meV. 
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Fig. 9. (color online) The basis of Wannier orbitals in FeSb2, presented as an isosurface plot in real space within the coordinate 
system of the primitive unit cell. In the first row (from left to right), we show one of the orbitals with predominant e g 
character (pointing approximately in the direction of the ligands) and two of the orbitals with mainly t2g character, respectively. 
In the second row we show, from left to right, the p x , Py, Pz orbitals. 




G ZT YGXS RU G ZT YGXS RU G ZT YGXS RU 




G ZT YGXS RU G ZT YGXS RU G ZT YGXS RU 



Fig. 10. Partial contributions to the bandstructure corresponding to the Wannier orbitals visualized in Figure 9. From top left 
to bottom right, the widths of the bands indicate the contribution of the e g , £2g(Sg ), i2g(aig), Px, Py, and p z Wannier function, 
respectively. Note the huge p-d hybridization in the energy region [-4 eV, 4 eV]. 
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4. Summary 



We have presented the implementation of an interface between the FLAPW code Wien2K and wannier90 
software for the construction of the maximally localized Wannier functions. The rational of this development 
is to provide a link between two widely used packages in the electronic structure community. We tried to 
provide examples of what we consider typical application of such a construction such as the construction of 
tight-binding Hamiltonians for complex systems, unfolding of complicated bandstructures or visualization. 
We also considered as an applications that have been so far less common such as strongly spin-orbit coupled 
Wannier orbitals. Last but not least we point out that the overlap matrices [Eq. (4)] can find application of 
its own, e.g. concerning the calculation of the interaction with an electro-magnetic field beyond the dipolc 
approximation. 
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